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Abstract. We present the first results from simulations of processes leading to planet formation in protoplanetary disks with 
different metallicities. For a given metallicity, we construct a two-dimensional grid of disk models with different initial masses 
and radii (Mo, Rq). For each disk, we follow the evolution of gas and solids from an early evolutionary stage, when all solids 
are in the form of small dust grains, to the stage when most solids have condensed into planetesimals. Then, based on the core 
accretion - gas capture scenario, we estimate the planet-bearing capability of the environment defined by the final planetesimal 
swarm and the still evolving gaseous component of the disk. We define the probability of planet-formation, P p , as the normalized 
fractional area in the (Mo, log^o) plane populated by disks that have formed planets inside 5 AU. With such a definition, and 
under the assumption that the population of planets discovered at R < 5 AU is not significantly contaminated by planets that 
have migrated from R > 5 AU, our results agree fairly well with the observed dependence between the probability that a star 
harbors a planet and the star's metal content. The agreement holds for the disk viscosity parameter a ranging from 10~ 3 to 10~ 2 , 
and it becomes much poorer when the redistribution of solids relative to the gas is not allowed for during the evolution of model 
disks. 

1. Introduction 

Well over a hundred extrasolar planets have already been cata- 
logued. It is very likely that radial velocity surveys performed 
in the solar neighbourhood have discovered a significant frac- 
tion of all Jupiter-like planets with periods shorter than ~10 yr 
(or, equivalently, with orbital radii smaller than ~5 AU), whose 
parent stars exhibit periodic radial velocity variations with am- 
plitudes > 10 m s _1 Such a collection of planets provides a 
good data set with which predictions of differ ent theo ries of 
planet formation can be compared (e.g. llda & Linl2004l) . 

Current observations strongly suggest that planet-bearing 
stars tend to have higher metallicities than field s tars (Fig , [p. 
This correlation was reported f or the first time b y Santos et alJ 
<200(J, and later confirmed bv lFischer & Valentil 12003I) . The 
effect finds a natural explanation in the core accretion - gas 
capture (CAGC) model for giant planet formation. The model 
predicts a slow accretion of planetesimals onto a protoplanetary 
core until the mass of the core reaches a few M ffi . Then, gas ac- 
cretion becomes significant, and around the core an extended 
envelope is formed, whose mass increases faster than that of 
the core. Eventually, the envelope becomes more massive than 
the core, and a runaway accretion of gas ensues, which termi- 
nates when tidal effects become important or the protoplane- 
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Fig. 1. The probability of planet occurrence as a function of the 
star's iron abundance (on the standard logarithmic scale with 
zero representing the solar value). Reproduced with permission 
from lFischer & Valentil J2003I) 



tary disk is dissipated. The length of the core accretion phase 
decrease s with increasing s urface density of the planetesimal 
swarm (Polla cket aljfl9 961. which, in turn, is an increasing 
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function of the original metal content of the protoplanetary disk 
and its central star. Consequently, the giant planets are expected 
to form more easily in disks with higher metallicities. 

However, dust evolves due to different physical 
processes than the gaseous component of the disk 
( Weidenschilling & Cuzzi 1993) As a result, a significant 
redistribution of solids takes place with respect to the gas, 
and, locally, the surface density of solids can be consider- 
ably e nhanced compared to the initial one (Weidenschilling 
20031 IStepinski & Valageaslll997l) . Thus, statistical tests of 
the CAGC scenario must be based on evolving models of 
protoplanetary disks rather than static configurations which 
only serve as mass supply to planetary cores and envelopes. 
Such models must include the global evolution of solid s. An 
extensive survey of such models was reported by Kornet et al. 
J200ll) . and applied to the pro blem of extrasolar planets 
formation by iKornet et alJ J2002I) . In the latter paper it was 
found that effects of dust evolution allow for in situ formation 
of giant planets as close as 2 AU from the star. Here we obtain 
a much larger sample of disk models with different primordial 
metallicities, estimate how many of them could form planets, 
and compare the estimates to the observed correlation between 
stellar iron abundance and the probability of planet occurrence. 

In §2 we briefly explain our approach to the evolution of a 
protoplanetary disk and planet formation. The results are pre- 
sented in §3, while in §4 we summarize them and critically 
assess the simplifications on which our method is based. 

2. Method of calculation 

2.1. The disk 

We model the protoplanetary disk as a two component fluid, 
consisting of gas and solids. T he gaseous com ponent is de- 
scribed by the analytical model o f lSter>inskildl998b which gives 
the surface density of gas, Y. g , as a function of radius r and 
time t, in terms of a solution to the viscous diffusion equation. 
The viscosity coefficient is given by the standard a prescrip- 
tion. All other quantities characterizing the gas are calculated 
from a thin disk approximation, assuming vertical thermal bal- 
ance. To take into account the effect of different metallicities 
on the structure of the gaseous disk, we scaled the opacity in 
Stepinski's model by a constant factor Z, independent of r or t 
and representing the metallicity of the disk expressed in solar 
units. Such a simple-minded approach is justified by the fact 
that the opacity in protoplanetary disks is mainly due to dust 
grains and molecules. 

The main assumptions underlying our approach to the evo- 
lution of solids are (1) at each radius the particles have the 
same size (which in general varies over time), (2) there is only 
one component of dust, in this case corresponding to high- 
temperature silicates with the evaporation temperature T evap = 
1350 K and bulk density 6.3 g cirT 3 (this choice is justified 
by the fact that in most cases the surface density of a water- 
ice planetesimal swarm would be too low in the range of dis- 
tances considered here, i.e. r < 5 AU, to enab le form ation 
of a g iant planet in a reasonable time; see also IKornet et alJ 
( 2004)), (3) all collisions between particles lead to coagulation, 



(4) when the temperature exceeds T evap , local solids immedi- 
ately sublimate and the vapor evolves at the same radial veloc- 
ity as the gas component, (5) when the disk temperature falls 
below T emp , local vapor immediately condenses into grains 
with a radius of 10~ 3 cm (minimum allowed by the numeri- 
cal code), (6) the radial velocity of solid particles is entirely 
determined by the effects of gas drag, (7) the relative velocities 
of solid particles at which they collide are computed accord- 
ing to the turbulent model described by Stepin ski & Valag eas 
( 1997), (8) the evolution of solids does not affect the evolu- 
tion of gas. At each radius, the vertical extent of the solid par- 
ticle distribution is calculated and is evolved in time, so the 
effect of sedimentation of solids toward the midplane of the 
disk is taken int o account. All assumptions and approximations 
are discussed in Stepinski & Valageas ( 1997), Stepinski ( 1998) 
and lKornet et al.N200ll) 

The evolution of solids is governed by two equations. The 
first of them is the continuity equation for surface density of 
solid material, The second one, describing the evolution of 
grain sizes, can be interpreted as the continuity equation for 
size-weighted surface density of solids, E fl = a(r)E„ where 
a(r) is the radius of solid particles at a distance r from the star. 
The equations are solved numerically on a moving grid whose 
outer edge follows the outer e dge of the dust dis k. The details 
of the method can be found in IKornet et al.l |2001 ). 

We parameterize the initial conditions by the initial mass 
of the disk, Mq (in M Q ) and the decimal logarithm of the initial 
outer radius of the disk, Rq (in AU). From these two quantities 
the initial total angular momentum of the disk jo, which we 
used in our previous calculations, can be calculated from the 
formula 



jo = 6.9M R l /2 M l J 2 



(D 



in units of 10 52 g cm 2 s _1 , where M* is the mass of t he ce ntral 
star expressed in M Q (see also eq. (211in lstepinskil dl998h . In 
the present survey M+ = 1 . Once Mo and Rq are specified, the 
initial distribution of gas surface density E g o = S g (r, t = 0), as 
well as the complete evoluti on of the gaseou s disk, are obtained 
from the analytical model of Stepinski ( 1993). 

We assume that at t — all solids are in the form of grains 
with radii a = 10~ 3 cm, well mixed with the gas. Thus, the 
ratio of 2 S o to S g o is independent of r, and proportional to disk 
metallicity, Z. For a given metallicity, the initial surface density 
distribution of solids is 



£ j0 (r) = 6x 10 _3 ZX 



•*.<>' 



(2) 



where Z is expressed in solar units. 
2.2. The planets 

We model the formation of a giant planet in situ, i.e. we do 
not include effects of migration, and the orbital parameters 
of the planet do not vary in time. Its evoluti on is followed 
acc ording to the approach dev eloped bv IPollack et alJ {1 996) 
and iBodenheimer et alJ J2000l) . The planet consists of a solid 
core and a gaseous envelope, both of which accrete mass from 
the disk. Our implementation of the evolutionary procedure is 
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based on the following assumptions: (1) core accretion starts 
when all solids are in the form of planetesimals with the same 
radius a — 2 km, (2) the surface density of solids in the feed- 
ing zone decreases with time as material accretes onto the core, 
(3) at each time the planetesimals are well mixed through the 
feeding zone of the planet; thus E s is always uniform in space, 
but usually decreasing with time, (4) the planetesimals do not 
migrate into the feeding zone from outside, or vice versa, but 
they can be overtaken by the boundary of the feeding zone as it 
expands due to the growing mass of the planet. 

With these assumptions, the minimum surface density of 
the planetesimal swarm 2 Jjmm , needed to form a Jupiter-mass 
planet in less than 3 Myr can be determined as a function of dis- 
ta nce r from the star. T he details of the procedure are discussed 
in K ornet et al.l d2002l) : here we only note that 2-km planetes- 
imals are on essentially Keplerian orbits, i.e. drag force has 
negligible effects on their distribution within the gaseous disk. 
The limit of 3 Myr is supported by recent ob servations (e.g. 
ISicilia-A guilar et aT1 l2004t iHaisch et alJEoQll) . and it leaves a 
reasonable safety margin for the planet to complete its evolu- 
tion before the canonical time of disk dispersal (10 7 yr) elapses. 
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Fig. 2. The rate of planet occurrence as a function of the pri- 
mordial metallicity of protoplanetary disks. Lines marked with 
different symbols are obtained for different values of a, as la- 
beled in the upper left corner. Filled symbols: models with the 
redistribution of solids; open symbols: models with a constant 
Ss/Sg ratio. The histogram shows the observational data. 



3. Results 

To cover the range of masses and sizes of protoplanetary disks 
indicated by observations, we consider 1 1 values of Mo uni- 
formly distributed between 0.02 and 0.2, and 10 values of Rq 
so chosen that lg Rq is uniformly distributed between 1 and 4. 
We evolve this set of models for seven different metallicities 
(Z = 0.2, 0.5, 0.75, 1, 1.5, 2.0, 2.5 and 3.0), and three values 
of the viscosity coefficient {a = 10~ 3 , 10~ 2 and 10~'). Using 
the procedure described in £12. II we follow each model until 
all solids are in the form of "Keplerian" planetesimals or they 
have accreted onto the star. Next, each model with planetesi- 
mals is scanned to determine whether anywhere in the disk the 
surface density of the final planetesimal swarm, exceeds 
£j,min( r )' Models with a range of r such that ~L s j(r) > 2 Jim ,„(r) 
are planet bearing; however for further analysis only those that 
form planets closer than 5 AU from the star are selected. We 
introduce this limitation in order to account for the fact that (as 
we already mentioned in §1) the sample of observed extrasolar 
planets is reasonably complete to ~5 AU only. 

For given a and Z our grid of Mo and lg Rq may be thought 
of as a nonrandom sample drawn from the unknown ensem- 
ble of disk masses and radii. Let "Rp t s be the region in the 
[Mo,lg^o] plane occupied by models that form planets at r < 
5 AU from the star, and let A p $ be the area of < R P: 5- The rate of 
planet occurrence within the sample is given by the ratio 



where A is the area of the region occupied by the whole grid of 
models. 

The results are presented in Fig. |2 It shows the rate of 
planet occurrence as a function of disk metallicity for different 
values of a. As one might expect, P p is an increasing function 
of Z: in more metal-rich disks there is more solid material, so 



the planets form more easily. We believe that deviations from 
that rule result from the poor resolution of our grid of models. 

There is also a monotonic trend seen when P p is consid- 
ered as a function of a at constant Z: giant planets form more 
easily in less viscous disks. This is because for every a and ev- 
ery Z the lower border of the planet-bearing region runs very 
close, and nearly parallel, to the line X that separates models in 
which all solids accrete onto the star from those in which some 
solids survive (see Fig. [3}- The observed trend results from two 
competing factors. On the one hand, to a good approximation 
£, marks the location of models in which the ini tial tempera- 
ture at the o uter edge of the disk is equal to T evap jKornet et alJ 
2001 , 2002). Among disks with the same Mq and Rq, those with 
smaller a are cooler as a result of lower energy generation rate. 
Thus, when a decreases, £, must move toward smaller values of 
Rq, and the area of the planet-bearing region increases. On the 
other hand, decreasing a reduces the maximum Rq for which 
planet formation is possible in disks with the same Mo. This is 
because in cooler disks the solids gain smaller inward velocities 
dWeidenschillingll977lEornet et all20 04). As a result, for the 
same disk parameters except for a lower a, the redistribution of 
solids is less efficient, and the final planetesimal swarm is more 
extended. A more extended swarm implies lower surface den- 
sity of solids and less favorable conditions for planet formation. 
In particular, when we move on the [Mo, \gRo\ plane along an 
Mo = const line such that Rq increases, we leave the planet- 
bearing region sooner (in other words, decreasing a causes the 
upper boundary of the planet-bearing region to move down- 
wards on the [Mq,Rq\ plane). Our calculations show that the 
first of these two effects is the dominant one. Therefore, for 
a given metallicity the fraction of disks with giant planets in- 
creases with decreasing a. 

To illustrate the importance of the redistribution of solids 
for the process of planet formation, we calculated P p for the 
same grid of initial conditions, but based on 2 ,,o instead of S s j. 
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Fig. 3. As functions of the initial disk mass Mq (in solar 
masses) and the logarithm of the initial radius Rq in astronom- 
ical units, the contours and grey scale, respectively, give the 
inner and outer radius of the region around the central star 
where giant planet formation is possible in a maximum of 3 
Myr. White circles indicate disk models where the solid sur- 
face density is everywhere below the critical value for planet 
formation. Black circles indicate disks in which all of the solid 
material accretes onto the star. The disk viscosity parameter is 
a = 1 x 10 -2 , and the primordial metallicity is solar. 

The results are shown in Fig.|5] As expected, the redistribution 
of solid material is most important for a = 10 _I , but even for 
a = 10" 3 it cannot be neglected. For example, with a = 10~ 3 
and Z = 3, including redistribution increases the ratio of planet 
occurrence by a factor as 1.5. For lower metallicities the effect 
is even stronger. In particular, it causes Z m ,„, the minimum Z 
at which planets can form, to decrease from- 0.7 to ~ 0.4 and 
from ~ 1.2 to ~ 0.85, respectively, for a = 10 -3 and 10 _1 . It can 
be seen that even in low-metallicity disks 2 S can locally grow 
above ,„,„, creating conditions favorable for planet formation. 

With a = 10~ 3 and the redistribution of solids included, the 
calculated rate of occurrence of planets is in good agreement 
with observational data. Obviously, it is a mere coincidence, 
without any practical meaning. Our rates strongly depend on 
the extent and position of the grid of initial conditions in the 
[Mo,Rq] plane. Suppose that we increase the maximum \gRo 
above the considered value of 4. The rate of planet occurrence 
would decrease, because we would add only models in which 
li s j is not high enough to satisfy the requirements for planet 
formation. To make our results independent of this effect, for 
every value of a we normalized P p {Z) in such a way as to repro- 
duce the observational value of ~ 0.2 for Z » 3. The normal- 
ized rates are shown in Fig. [4] Now, curves for both a = 10~ 3 
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Fig. 4. Rate of planet occurrence as a function of the primordial 
metallicity in models with the redistribution of solids included. 
Same as in Fig. [2] but all curves are normalized to match the 
observational data for lg Z = 0.4. 
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Fig. 5. Rate of planet occurrence as a function of the primor- 
dial metallicity in models without the redistribution of solids. 
Same as in Fig. [2] but all curves are normalized to match the 
observational data for lg Z = 0.4. 

and 10~ 2 fall in the range of observational values, while disks 
with a - 10 _1 still produce too few planets for small values of 
Z. For comparison, in Fig.[5]we shows the normalized rates ob- 
tained from models without the redistribution of solids. While 
for a = 10~ 3 and 10~ 2 they are consistent with the observa- 
tional data for metal rich disks, they are significantly lower 
for small values of Z, and they give too high values of Z mi „. 
For a = 10~', the results do not match observations except for 
Z ~ 3, where the agreement is enforced by the normalization 
procedure. 

4. Discussion and conclusions 

Based on a simple approach to the evolution of solids in proto- 
planetary disks we calculated rates of giant planet occurrence 
on orbits smaller than 5 AU. We find that in more metal-rich 
disks the planets form more easily, and we qualitatively repro- 



K. Kornet P. Bodenheimer, M. Rozyczka, T.F. Stepinski: Giant planets in disks with different metallicities 



5 



0.2 



0.15 - 



0.1 



0.05 



observations 
oe=10-I 
a=10"^ 
a=10 



-0.8 -0.6 -0.4 -0.2 

igz 



0.2 0.4 



Fig. 6. Same as Fig.|4] but based on disk models with initial 
parameters uniformly distributed on the [Mo, lg jo] plane. 



duce the observed correlation between the probability that a gi- 
ant planet would be found at a star and the metal content of that 
star. We are able to match this correlation quantitatively for disk 
models with the viscosity parameter a ranging from ~ 10~ 3 to 
~ 10~ 2 . We also show that in the low-metallicity regime the 
calculated occurrence rates agree with the observed ones only 
if the redistribution of solids is included in the models. In such 
models the surface density of solids locked in the final plan- 
etesimal swarm can be locally much higher than the initial sur- 
face density of the dust. This effect is particularly important if 
we want to match the observed value of limiting stellar metal- 
licity Z m ,„ below which no planets have been detected. Again, 
only models which allow for the redistribution of solids are suc- 
cessful, while models which do not include this effect produce 
Z m i„ which is significantly too high. This finding supports our 
approach to the evolution of protoplanetary disks, whose key 
ingredient is the radial drift of solids due to the gas drag. 

The above conclusions are based on the assumption that 
our grid of models is a representative sample of the real distri- 
bution P(Mo,/?o) of initial masses and radii of protoplanetary 
disks, which is not known. Such an assumption is acceptable 
provided that the real disks are uniformly (or almost uniformly) 
distributed on the [Mo,lgRo] plane. Adopting it we demand 
that very large disks be relatively rare, which seems to be an 
entirely reasonable requirement (note that P(Mq,Ro) must ap- 
proach for Ro — > oo). We tested the sensitivity of our results 
to the way the initial parameters are sampled by repeating the 
calculations with initial models uniformly distributed on the 
[Mo,lgjo] plane. The results, normalized to the observational 
data the same way as before, are shown in Fig. [6] The agree- 
ment is poorer, but the same general trend is still well visible. 
Moreover, we see that Z,„,„ does not depend on the sampling 
procedure (this conclusion is restricted to samples which do 
not include extremely massive disks with Mq > 0.2 M Q ; note 
however that such disks would most likely be gravitationally 
unstable). 

Our results are based on a highly simplified scenario of 
the evolution of solids in protoplanetary disks. The assump- 
tion concerning the size distribution of solid particles, ac- 



cording to which at each distance from the star all parti- 
cles have the same diameter, se ems to be particularly radi- 
cal. H owev er, detailed work bv | Morfilll (|1985) Mi zuno et alJ 
lll988t) and IWeidenschilling l lll997h showed that the distribu- 
tion of solids quickly converges to a stage in which most of the 
mass is concentrated in a narrow range of sizes at the momen- 
tary maximum size. Their conclusion has recently been corrob- 
orated by the behaviour of solids in the two-dimensional disk 
model by Weidenschilling (2003). Thus, despite its simplicity, 
our our approach to the size disitribution function can be re- 
garded as a reasonable approximation. 

Our next radical assumption is the 100% efficiency of co- 
agulation In reality, one hardly expects that collisions between 
solid particles always lead to stickin g without fragmen tation. 
However, calculations reported by iKornet et al.l (12.0041) show 
that if the main process responsible for particle growth is colli- 
sional coagulation, then its efficiency must be high - otherwise 
the disk would rapidly lose most of its solid material, disabling 
planet formation or severely diminish ing its possibility. Note 
also that recent observations by Chakrabor tv et alj j2004) indi- 
cate that asteroid-sized bodies can exist in disks as young as 0. 1 
Myr. Again, if they have grown due to coagulation, the sticking 
efficiency must have been very high. 

Migration of planets due to gravitational interactions with 
the disk is not included in our approach. Obviously, we do not 
claim that the process of migration does not operate. On the 
contrary, we have to assume some form of redistribution of 
planets in order to explain the origin of "hot Jupiters", whose 
orbits are much tighter than the tightest orbits allowed by our 
models (~ 1 AU). However, as long as we are interested only 
in global rates of planet occurrence, and not in the distribution 
of orbital parameters, our results remain valid provided that the 
number of planets which migrate to r < 5 AU from larger or- 
bits is negligible compared to the number of planets which are 
born there. 

The fact that such a simple model can reproduce the avail- 
able data has two possible interpretations. According to the 
pessimistic one we observe a mere coincidence which does not 
have any physical meaning. According to the optimistic one 
we identify the most important process governing the evolution 
of solids in a protoplanetary disk (radial drift, induced by the 
gas drag, and associated with collisional coagulation), and we 
provide support for the core accretion - envelope capture sce- 
nario of giant planet formation. More sophisticated disk models 
based on fewer assumptions, which we are presently working 
on, should help us decide which of the two is true. 
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